SARS-CoV-2 spike opening dynamics and energetics reveal the individual roles of glycans and their collective impact

The trimeric spike (S) glycoprotein, which protrudes from the SARS-CoV-2 viral envelope, binds to human ACE2, initiated by at least one protomer’s receptor binding domain (RBD) switching from a "down” (closed) to an "up” (open) state. Here, we used large-scale molecular dynamics simulations and two-dimensional replica exchange umbrella sampling calculations with more than a thousand windows and an aggregate total of 160 μs of simulation to investigate this transition with and without glycans. We find that the glycosylated spike has a higher barrier to opening and also energetically favors the down state over the up state. Analysis of the S-protein opening pathway reveals that glycans at N165 and N122 interfere with hydrogen bonds between the RBD and the N-terminal domain in the up state, while glycans at N165 and N343 can stabilize both the down and up states. Finally, we estimate how epitope exposure for several known antibodies changes along the opening path. We find that the BD-368-2 antibody’s epitope is continuously exposed, explaining its high efficacy.

where p(λ ,t|λ 0 , 0) is the probability of finding the system at λ after time t, given it was at λ 0 at time 0, D(λ ) is a position-dependent diffusion constant, and F(λ ) is the free energy at λ . By rearranging the terms, the MFPTτ FP , or the rate inverse k −1 , from the initial (A) to the final (B) state is given by: While F(λ ) was readily available from the PMF obtained through REUS, D(λ ) was approximated using a generalized-Langevin-equation-based method, derived by Roux and co-workers and implemented by Gaalswyk et al. 6 . From a time series of a coordinate x k with the system simulated under a harmonic restraint, the method computes the diffusion constant D k along x k by relating it to its velocity autocorrelation function (VACF). With a series of coordinate transformations 7 , D k is transformed from the Cartesian space to the collective variable space (D i j ), and then to the path variable space (D(λ )), which was then inserted back into Eq. S2.
To determine D k for each atom and each of its coordinates involved in the definition of the collective variables d and ϕ , we ran a 1-ns simulation for each window along the MEP, with all C α atoms of the protein restrained by a force constant of 5 kcal mol −1 Å −2 . A 2-fs time step was used without the application of HMR. Other simulation parameters were identical to the REUS simulations.
Kinetics analysis. We modeled the down-to-up transition and subsequent binding of the RBD to ACE2 according to the chemical equation The associated master equation is The rates k open and k close come from our own calculations. Hydrogen bond analysis. The number of hydrogen bonds formed between RBD-A and other domains of the spike in the REUS trajectories were measured using the HBonds Plugin of VMD. When the distance between a hydrogen-bond donor atom (D) and an acceptor atom (A) is below 3.5 Å and the angle D-H-A is less than 35 • from 180 • , a hydrogen bond was considered to be formed between D and A. The d-ϕ space is broken into small bins and the number of hydrogen bonds was averaged over each bin.
Contact analysis. Contact calculations were performed using a cutoff distance of 3.5 Å. In other words, when two heavy atoms from two different selections come within 3.5 Å, we count that as one contact. For the analysis done in Figs. 4c-d and S8, the d-ϕ space is broken into small bins and the number of contacts was averaged over each bin. For the analysis done in Fig. S9, we separated conformations from the MEP in two categories and the number of contacts were averaged separately for each category. The two categories were defined by the d values.
The conformations with d ≤ 54.9 Å were considered as the down state, while the rest were taken as the up state.
Antibody accessible surface area (AbASA). AbASA calculations were performed using the solvent-accessible surface area measurement tool as implemented in VMD with a 7-Å probe for each frame along the MEP. For the glycosylated system, we performed two accessible surface area calculations for each frame around each epitope using the -restrict option: AbASA of (calculation #1) protein and (calculation #2) protein + glycans. The difference in the accessible surface area obtained from these two calculations gives the coverage provided by the glycans only.